## REPLICATION CODE FOR: #####################################################################################
#
# Guess, Nagler, and Tucker. "Less than you think: Prevalence and predictors of fake news dissemination
#                             on Facebook." Science Advances 5, eaau4586 (2019).
# DOI: 10.1126/sciadv.aau4586
#
##############################################################################################################

library(tidyverse)

setwd("")
load("gnt_replication_data.RData")

## TABLE 1

library(xtable)

xtable(table(cut(dat$num_fake_shares, breaks = c(-1, 0, 1, 2, 3, 4, 10, max(dat$num_fake_shares, na.rm = TRUE)), labels = c("0", "1", "2", "3", "4", "5-10", "11-50"))))
round(table(cut(dat$num_fake_shares, breaks = c(-1, 0, 1, 2, 3, 4, 10, max(dat$num_fake_shares, na.rm = TRUE)), labels = c("0", "1", "2", "3", "4", "5-10", "11-50")))/1191, digits = 3)

## TABLE 2

model.cp1a <- glm(num_fake_shares ~ ideology + age + female + black + educ + faminc, 
                  weights = weight_svypop_w3, data = dat, family = "quasipoisson")

model.cp1 <- glm(num_fake_shares ~ ideology + age + female + black + educ + faminc + num_posts, 
                 weights = weight_svypop_w3, data = dat, family = "quasipoisson")

model.cp2a <- glm(num_fake_shares_ag ~ ideology + age + female + black + educ + faminc, 
                  weights = weight_svypop_w3, data = dat, family = "quasipoisson")

model.cp2 <- glm(num_fake_shares_ag ~ ideology + age + female + black + educ + faminc + num_posts, 
                 weights = weight_svypop_w3, data = dat, family = "quasipoisson")

library(stargazer)

varnames <- c("(Intercept)", "Very Liberal", "Liberal",
              "Moderate", "Conservative", "Very Conservative", 
              "Age: 30-44", "Age: 45-65", "Age: Over 65", "Female", "Black", 
              "Education", "Income", "Number of links shared") 

stargazer(model.cp1a, model.cp1, model.cp2a, model.cp2,
          style = "apsr", align = TRUE, 
          title = "Determinants of Fake News Sharing on Facebook",
          covariate.labels = c(varnames[-1]),
          dep.var.labels = c("Number of Stories Shared", "Number of Stories Shared (A\\&G)"), 
          keep.stat = c("n", "null.dev", "res.dev"), 
          notes = "Quasi-Poisson models.")

## TABLE 3

# Column 1 (all):

round(prop.table(table(dat$pid3, useNA = "always")), digits = 2) # 1 = Dem
mean(dat$ideo5, na.rm = TRUE)
round(prop.table(table(dat$hrc, useNA = "always")), digits = 2)
mean(dat$voted, na.rm = TRUE)
mean(dat$knowl, na.rm = TRUE)
mean(dat$ageyrs, na.rm = TRUE)
round(prop.table(table(dat$educ, useNA = "always")), digits = 2) # 1 and 2 = HS or less


# Column 2 (report FB):

round(prop.table(table(dat$pid3[which(dat$fb_yes == 1)], useNA = "always")), digits = 2)
mean(dat$ideo5[which(dat$fb_yes == 1)], na.rm = TRUE)
round(prop.table(table(dat$hrc[which(dat$fb_yes == 1)], useNA = "always")), digits = 2)
mean(dat$voted[which(dat$fb_yes == 1)], na.rm = TRUE)
mean(dat$knowl[which(dat$fb_yes == 1)], na.rm = TRUE)
mean(dat$ageyrs[which(dat$fb_yes == 1)], na.rm = TRUE)
round(prop.table(table(dat$educ[which(dat$fb_yes == 1)], useNA = "always")), digits = 2)

round(prop.table(table(dat$fb_post[which(dat$fb_yes == 1)], useNA = "always")), digits = 2) # 6 = Several times a day
round(prop.table(table(dat$lookoften[which(dat$fb_yes == 1)], useNA = "always")), digits = 2)

table(dat$fb_yes) # N


# Column 3 (provided FB access):

round(prop.table(table(dat$pid3[which(dat$fb_data_yes == 1)], useNA = "always")), digits = 2)
mean(dat$ideo5[which(dat$fb_data_yes == 1)], na.rm = TRUE)
round(prop.table(table(dat$hrc[which(dat$fb_data_yes == 1)], useNA = "always")), digits = 2)
mean(dat$voted[which(dat$fb_data_yes == 1)], na.rm = TRUE)
mean(dat$knowl[which(dat$fb_data_yes == 1)], na.rm = TRUE)
mean(dat$ageyrs[which(dat$fb_data_yes == 1)], na.rm = TRUE)
round(prop.table(table(dat$educ[which(dat$fb_data_yes == 1)], useNA = "always")), digits = 2)

round(prop.table(table(dat$fb_post[which(dat$fb_data_yes == 1)], useNA = "always")), digits = 2)
round(prop.table(table(dat$lookoften[which(dat$fb_data_yes == 1)], useNA = "always")), digits = 2)

table(dat$fb_data_yes)


# Column 4 (FB data linked to survey):

round(prop.table(table(dat$pid3[which(!is.na(dat$num_fake_shares))], useNA = "always")), digits = 2)
mean(dat$ideo5[which(!is.na(dat$num_fake_shares))], na.rm = TRUE)
round(prop.table(table(dat$hrc[which(!is.na(dat$num_fake_shares))], useNA = "always")), digits = 2)
mean(dat$voted[which(!is.na(dat$num_fake_shares))], na.rm = TRUE)
mean(dat$knowl[which(!is.na(dat$num_fake_shares))], na.rm = TRUE)
mean(dat$ageyrs[which(!is.na(dat$num_fake_shares))], na.rm = TRUE)
round(prop.table(table(dat$educ[which(!is.na(dat$num_fake_shares))], useNA = "always")), digits = 2)

round(prop.table(table(dat$fb_post[which(!is.na(dat$num_fake_shares))], useNA = "always")), digits = 2)
round(prop.table(table(dat$lookoften[which(!is.na(dat$num_fake_shares))], useNA = "always")), digits = 2)

sum(table(dat$num_fake_shares))


# T-tests: report vs. provide access

dat$fb_post2 <- ifelse(dat$fb_post == 6, 1, 0)

t.test(dat$dem[which(dat$fb_yes == 1)], dat$dem[which(dat$fb_data_yes == 1)])
t.test(dat$ideo5[which(dat$fb_yes == 1)], dat$ideo5[which(dat$fb_data_yes == 1)])
t.test(dat$hrc[which(dat$fb_yes == 1)], dat$hrc[which(dat$fb_data_yes == 1)])
t.test(dat$voted[which(dat$fb_yes == 1)], dat$voted[which(dat$fb_data_yes == 1)])
t.test(dat$knowl[which(dat$fb_yes == 1)], dat$knowl[which(dat$fb_data_yes == 1)])
t.test(dat$ageyrs[which(dat$fb_yes == 1)], dat$ageyrs[which(dat$fb_data_yes == 1)])
t.test(dat$hs[which(dat$fb_yes == 1)], dat$hs[which(dat$fb_data_yes == 1)])

t.test(dat$fb_post2[which(dat$fb_yes == 1)], dat$fb_post2[which(dat$fb_data_yes == 1)])
t.test(dat$lookoften[which(dat$fb_yes == 1)], dat$lookoften[which(dat$fb_data_yes == 1)])


## FIGURE 1

library(hrbrthemes)

# Left panel

p <- ggplot(filter(dat, party %in% c("Democrat", "Republican", "Independent")), aes(num_posts)) + geom_histogram(bins = 50)
pq <- p + ggtitle("Distribution of Fake News Shares") + ylab("Number of respondents") + 
  xlab("Number of web links shared") + theme(legend.position = "bottom", legend.title = element_blank()) + labs(title = NULL)

pq + theme_ipsum_rc()

# Right panel

p <- ggplot(filter(dat, party %in% c("Democrat", "Republican", "Independent")), aes(num_fake_shares, fill = party)) + geom_histogram(bins = 50)
pq <- p + ggtitle("Distribution of Fake News Shares") + ylab("Number of respondents") + 
  scale_fill_manual("", values = c("red", "blue", "darkgray")) +
  xlab("Number of fake news stories shared") + theme(legend.position = "bottom", legend.title = element_blank()) + labs(title = NULL)

pq + theme_ipsum_rc()

## FIGURE 2

library(spatstat)

# (a)

dat$party <- relevel(dat$party, ref = "Other/Not sure")
dat$party <- relevel(dat$party, ref = "Independent")
dat$party <- relevel(dat$party, ref = "Democrat")
dat$party <- relevel(dat$party, ref = "Republican")
ggdat <- dat %>% filter(!is.na(party)) %>% group_by(party) %>% summarize(mean_fakes = mean(num_fake_shares, na.rm = TRUE),
                                                                         mean_fakes_w = weighted.mean(num_fake_shares, w = weight_svypop_w3, na.rm = TRUE),
                                                                         sd_fakes = sd(num_fake_shares, na.rm = TRUE)/sqrt(n()),
                                                                         sd_fakes_w = sqrt(weighted.var(num_fake_shares, w = weight_svypop_w3, na.rm = TRUE)/n()))

p <- ggplot(filter(ggdat, !is.na(party)), aes(party, mean_fakes_w))
pq <- p + geom_bar(stat = "identity") + ylim(-.02, 1.23) + geom_errorbar(aes(ymin = mean_fakes_w-1.96*sd_fakes_w, ymax = mean_fakes_w+1.96*sd_fakes_w), width = .1, col = "darkgray") +
  ggtitle("Party ID and Fake News Shares") + ylab("Mean number of fake news stories shared") + 
  xlab("Party identification") + theme(legend.position = "none") + labs(title = NULL)

pq + theme_ipsum_rc(grid = "Y")

# (b)

ggdat <- dat %>% filter(!is.na(age)) %>% group_by(age) %>% summarize(mean_fakes = mean(num_fake_shares, na.rm = TRUE),
                                                                     mean_fakes_w = weighted.mean(num_fake_shares, w = weight_svypop_w3, na.rm = TRUE),
                                                                     sd_fakes = sd(num_fake_shares, na.rm = TRUE)/sqrt(n()),
                                                                     sd_fakes_w = sqrt(weighted.var(num_fake_shares, w = weight_svypop_w3, na.rm = TRUE)/n()))

p <- ggplot(ggdat, aes(age, mean_fakes_w))
pq <- p + geom_bar(stat = "identity") + ylim(-.003, 1.23) + geom_errorbar(aes(ymin = mean_fakes_w-1.96*sd_fakes_w, ymax = mean_fakes_w+1.96*sd_fakes_w), width = .1, col = "darkgray") +
  ggtitle("Age and Fake News Shares") + ylab("Mean number of fake news stories shared") + 
  xlab("Age group") + theme(legend.position = "none") + labs(title = NULL)

pq + theme_ipsum_rc(grid = "Y")

# (c)

ggdat <- dat %>% filter(!is.na(ideology)) %>% group_by(ideology) %>% summarize(mean_fakes = mean(num_fake_shares, na.rm = TRUE),
                                                                               mean_fakes_w = weighted.mean(num_fake_shares, w = weight_svypop_w3, na.rm = TRUE),
                                                                               sd_fakes = sd(num_fake_shares, na.rm = TRUE)/sqrt(n()),
                                                                               sd_fakes_w = sqrt(weighted.var(num_fake_shares, w = weight_svypop_w3, na.rm = TRUE)/n()))

p <- ggplot(filter(ggdat, !is.na(ideology), ideology != "Not sure"), aes(ideology, mean_fakes_w))
pq <- p + geom_bar(stat = "identity") + ylim(-.02, 1.29) + geom_errorbar(aes(ymin = mean_fakes_w-1.96*sd_fakes_w, ymax = mean_fakes_w+1.96*sd_fakes_w), width = .1, col = "darkgray") +
  ggtitle("Ideology and Fake News Shares") + ylab("Mean number of fake news stories shared") + 
  xlab("Ideological self-placement") + labs(title = NULL)

pq + theme_ipsum_rc(grid = "Y")

# (d)

ggdat <- dat %>% group_by(posts) %>% summarize(mean_fakes = mean(num_fake_shares, na.rm = TRUE),
                                               mean_fakes_w = weighted.mean(num_fake_shares, w = weight_svypop_w3, na.rm = TRUE),
                                               sd_fakes = sd(num_fake_shares, na.rm = TRUE)/sqrt(n()),
                                               sd_fakes_w = sqrt(weighted.var(num_fake_shares, w = weight_svypop_w3, na.rm = TRUE)/n()))

p <- ggplot(filter(ggdat, !is.na(posts)), aes(posts, mean_fakes))
pq <- p + geom_bar(stat = "identity") + ylim(-.07, 1.2) + geom_errorbar(aes(ymin = mean_fakes-1.96*sd_fakes_w, ymax = mean_fakes+1.96*sd_fakes_w), width = .1, col = "darkgray") +
  ggtitle("Post Frequency and Fake News Shares") + ylab("Mean number of fake news stories shared") + 
  xlab("Number of links posted to Facebook") + labs(title = NULL)

pq + theme_ipsum_rc(grid = "Y")



# Supplementary Materials -------------------------------------------------


## TABLE S1

library(MASS)

model.nb1 <- glm.nb(num_fake_shares ~ ideology + age + female + black + educ + faminc + num_posts, 
                    weights = weight_svypop_w3, data = dat)
model.nb2 <- glm.nb(num_fake_shares_ag ~ ideology + age + female + black + educ + faminc + num_posts, 
                    weights = weight_svypop_w3, data = dat)

stargazer(model.nb1, model.nb2,
          style = "apsr", align = TRUE, 
          title = "Determinants of Fake News Sharing on Facebook",
          covariate.labels = c(varnames[-1]),
          dep.var.labels = c("Number of Stories Shared", "Number of Stories Shared (A\\&G)"), 
          omit.stat = c("f", "ser", "rsq"), notes = "Negative binomial regressions.")

## TABLE S2

model.ols1 <- lm(num_fake_shares ~ ideology + age + female + black + educ + faminc + num_posts, 
                 weights = weight_svypop_w3, data = dat)
model.ols2 <- lm(num_fake_shares_ag ~ ideology + age + female + black + educ + faminc + num_posts, 
                 weights = weight_svypop_w3, data = dat)

stargazer(model.ols1, model.ols2,
          style = "apsr", align = TRUE, 
          title = "Determinants of Fake News Sharing on Facebook",
          covariate.labels = c(varnames[-1]),
          dep.var.labels = c("Number of Stories Shared", "Number of Stories Shared (A\\&G)"), 
          omit.stat = c("f", "ser", "rsq"), notes = "OLS regressions.")

## TABLE S3

varnames <- c("(Intercept)", "Very Liberal", "Liberal",
              "Moderate", "Conservative", "Very Conservative", "Democrat",  "Republican", "Independent", "Other", 
              "Age: 30-44", "Age: 45-65", "Age: Over 65", "Female", "Black", 
              "Education", "Income", "Number of links shared")

dat$party <- factor(as.numeric(dat$pid3), labels = names(attr(dat$pid3, "labels"))[1:5])
dat$party <- relevel(dat$party, ref = "Not sure")

model.pid1a <- glm(num_fake_shares ~ ideology + party + age + female + black + educ + faminc + num_posts, 
                   weights = weight_svypop_w3, data = dat, family = "quasipoisson")
model.pid1b <- glm.nb(num_fake_shares ~ ideology + party + age + female + black + educ + faminc + num_posts, 
                      weights = weight_svypop_w3, data = dat)
model.pid1c <- lm(num_fake_shares ~ ideology + party + age + female + black + educ + faminc + num_posts, 
                  weights = weight_svypop_w3, data = dat)

stargazer(model.pid1a, model.pid1b, model.pid1c,
          style = "apsr", align = TRUE, model.names = FALSE,
          title = "Determinants of Fake News Sharing on Facebook",
          covariate.labels = c(varnames[-1]),
          dep.var.labels = "Number of Stories Shared", 
          omit.stat = c("f", "ser", "rsq"), notes = "Models: Quasi-Poisson, negative binomial, and OLS respectively.")

## TABLE S4

model.pid2a <- glm(num_fake_shares_ag ~ ideology + party + age + female + black + educ + faminc + num_posts, 
                   weights = weight_svypop_w3, data = dat, family = "quasipoisson")
model.pid2b <- glm.nb(num_fake_shares_ag ~ ideology + party + age + female + black + educ + faminc + num_posts, 
                      weights = weight_svypop_w3, data = dat)
model.pid2c <- lm(num_fake_shares_ag ~ ideology + party + age + female + black + educ + faminc + num_posts, 
                  weights = weight_svypop_w3, data = dat)

stargazer(model.pid2a, model.pid2b, model.pid2c,
          style = "apsr", align = TRUE, model.names = FALSE,
          title = "Determinants of Fake News Sharing on Facebook",
          covariate.labels = c(varnames[-1]),
          dep.var.labels = "Number of Stories Shared (A\\&G)", 
          omit.stat = c("f", "ser", "rsq"), notes = "Models: Quasi-Poisson, negative binomial, and OLS respectively.")

## TABLE S5

varnames <- c("(Intercept)", "Democrat",  "Republican", "Independent", "Other", 
              "Age: 30-44", "Age: 45-65", "Age: Over 65", "Female", "Black", 
              "Education", "Income", "Number of links shared") 

model.pid1a <- glm(num_fake_shares ~ party + age + female + black + educ + faminc + num_posts, 
                   weights = weight_svypop_w3, data = dat, family = "quasipoisson")
model.pid1b <- glm.nb(num_fake_shares ~ party + age + female + black + educ + faminc + num_posts, 
                      weights = weight_svypop_w3, data = dat)
model.pid1c <- lm(num_fake_shares ~ party + age + female + black + educ + faminc + num_posts, 
                  weights = weight_svypop_w3, data = dat)

stargazer(model.pid1a, model.pid1b, model.pid1c,
          style = "apsr", align = TRUE, model.names = FALSE,
          title = "Determinants of Fake News Sharing on Facebook",
          covariate.labels = c(varnames[-1]),
          dep.var.labels = "Number of Stories Shared", 
          omit.stat = c("f", "ser", "rsq"), notes = "Models: Quasi-Poisson, negative binomial, and OLS respectively.")

## TABLE S6

model.pid2a <- glm(num_fake_shares_ag ~ party + age + female + black + educ + faminc + num_posts, 
                   weights = weight_svypop_w3, data = dat, family = "quasipoisson")
model.pid2b <- glm.nb(num_fake_shares_ag ~ party + age + female + black + educ + faminc + num_posts, 
                      weights = weight_svypop_w3, data = dat)
model.pid2c <- lm(num_fake_shares_ag ~ party + age + female + black + educ + faminc + num_posts, 
                  weights = weight_svypop_w3, data = dat)

stargazer(model.pid2a, model.pid2b, model.pid2c,
          style = "apsr", align = TRUE, model.names = FALSE,
          title = "Determinants of Fake News Sharing on Facebook",
          covariate.labels = c(varnames[-1]),
          dep.var.labels = "Number of Stories Shared (A\\&G)",
          omit.stat = c("f", "ser", "rsq"), notes = "Models: Quasi-Poisson, negative binomial, and OLS respectively.")

## TABLES S7 + S8

dat$any_fake_shares <- (dat$num_fake_shares > 0)*1

mod.bin.1 <- lm(any_fake_shares ~ ideology + party + age + female + black + educ + faminc + num_posts, 
                weights = weight_svypop_w3, data = dat)
mod.bin.2 <- lm(any_fake_shares ~ party + age + female + black + educ + faminc + num_posts, 
                weights = weight_svypop_w3, data = dat)
mod.bin.3 <- lm(any_fake_shares ~ ideology + age + female + black + educ + faminc + num_posts, 
                weights = weight_svypop_w3, data = dat)

mod.bin.l1 <- glm(any_fake_shares ~ ideology + party + age + female + black + educ + faminc + num_posts, 
                  weights = weight_svypop_w3, data = dat, family = binomial(link = "logit"))
mod.bin.l2 <- glm(any_fake_shares ~ party + age + female + black + educ + faminc + num_posts, 
                  weights = weight_svypop_w3, data = dat, family = binomial(link = "logit"))
mod.bin.l3 <- glm(any_fake_shares ~ ideology + age + female + black + educ + faminc + num_posts, 
                  weights = weight_svypop_w3, data = dat, family = binomial(link = "logit"))


varnames <- c("(Intercept)", "Very Liberal", "Liberal",
              "Moderate", "Conservative", "Very Conservative", "Democrat",  "Republican", "Independent", "Other", 
              "Age: 30-44", "Age: 45-65", "Age: Over 65", "Female", "Black", 
              "Education", "Income", "Number of FB posts") 

stargazer(mod.bin.1, mod.bin.2, mod.bin.3, 
          style = "apsr", align = TRUE, model.names = FALSE,
          title = "Determinants of Fake News Sharing on Facebook",
          covariate.labels = c(varnames[-1]),
          dep.var.labels = "Shared Any Fake News Stories", 
          omit.stat = c("f", "ser", "rsq"), notes = "OLS models.")

stargazer(mod.bin.l1, mod.bin.l2, mod.bin.l3, 
          style = "apsr", align = TRUE, model.names = FALSE,
          title = "Determinants of Fake News Sharing on Facebook",
          covariate.labels = c(varnames[-1]),
          dep.var.labels = "Shared Any Fake News Stories", 
          omit.stat = c("f", "ser", "rsq"), notes = "Logistic regression models.")

## TABLE S9

model.cp1u <- lm(num_fake_shares_ag_url ~ ideology + age + female + black + educ + faminc + num_posts, 
                 weights = weight_svypop_w3, data = dat)
model.cp2u <- lm(num_fake_shares_ag_url_djt ~ ideology + age + female + black + educ + faminc + num_posts, 
                 weights = weight_svypop_w3, data = dat) # Trump
model.cp3u <- lm(num_fake_shares_ag_url_hrc ~ ideology + age + female + black + educ + faminc + num_posts, 
                 weights = weight_svypop_w3, data = dat) # Clinton

varnames <- c("Very Liberal", "Liberal",
              "Moderate", "Conservative", "Very Conservative", 
              "Age: 30-44", "Age: 45-65", "Age: Over 65", "Female", "Black", 
              "Education", "Income", "Number of links posted", "Constant")

stargazer(model.cp1u, model.cp2u, model.cp3u,
          style = "apsr", align = TRUE, 
          title = "Determinants of Fake News Sharing on Facebook",
          covariate.labels = c(varnames),
          dep.var.labels = c("Number of Stories Shared (A\\&G)", "Number of Stories Shared (A\\&G, pro-T)", "Number of Stories Shared (A\\&G, pro-C)"), 
          keep.stat = c("n", "null.dev", "res.dev"), 
          notes = "OLS models.")

## TABLE S10

model.bf1 <- glm(num_fake_shares_bf1 ~ ideology + age + female + black + educ + faminc + num_posts, 
                 weights = weight_svypop_w3, data = dat, family = "quasipoisson")
model.bf2 <- glm(num_fake_shares_bf1 ~ ideology + age + female + black + educ + faminc, 
                 weights = weight_svypop_w3, data = dat, family = "quasipoisson")
model.bf3 <- glm(num_fake_shares_bf1 ~ party + age + female + black + educ + faminc + num_posts, 
                 weights = weight_svypop_w3, data = dat, family = "quasipoisson")
model.bf4 <- glm(num_fake_shares_bf1 ~ party + age + female + black + educ + faminc, 
                 weights = weight_svypop_w3, data = dat, family = "quasipoisson")

varnames <- c("(Intercept)", "Very Liberal", "Liberal",
              "Moderate", "Conservative", "Very Conservative", "Democrat",  "Republican", "Independent", "Other", 
              "Age: 30-44", "Age: 45-65", "Age: Over 65", "Female", "Black", 
              "Education", "Income", "Number of FB posts") 

stargazer(model.bf1, model.bf2, model.bf3, model.bf4,
          style = "apsr", align = TRUE, 
          title = "Determinants of Fake News Sharing on Facebook",
          covariate.labels = c(varnames[-1]),
          dep.var.labels = c("Number of Stories Shared"), 
          omit.stat = c("f", "ser", "rsq"), notes = "Quasi-Poisson models.")

## TABLE S11

model.bf1 <- glm(num_fake_shares_bf2 ~ ideology + age + female + black + educ + faminc + num_posts, 
                 weights = weight_svypop_w3, data = dat, family = "poisson")
model.bf2 <- glm(num_fake_shares_bf2 ~ ideology + age + female + black + educ + faminc, 
                 weights = weight_svypop_w3, data = dat, family = "poisson")
model.bf3 <- glm(num_fake_shares_bf2 ~ party + age + female + black + educ + faminc + num_posts, 
                 weights = weight_svypop_w3, data = dat, family = "poisson")
model.bf4 <- glm(num_fake_shares_bf2 ~ party + age + female + black + educ + faminc, 
                 weights = weight_svypop_w3, data = dat, family = "poisson")

stargazer(model.bf1, model.bf2, model.bf3, model.bf4,
          style = "apsr", align = TRUE, 
          title = "Determinants of Fake News Sharing on Facebook",
          covariate.labels = c(varnames[-1]),
          dep.var.labels = c("Number of Stories Shared"), 
          omit.stat = c("f", "ser", "rsq"), notes = "Poisson models.")

## TABLE S12

model.os1 <- glm(num_fake_shares_os ~ ideology + age + female + black + educ + faminc + num_posts, 
                 weights = weight_svypop_w3, data = dat, family = "quasipoisson")
model.os2 <- glm(num_fake_shares_os ~ ideology + age + female + black + educ + faminc, 
                 weights = weight_svypop_w3, data = dat, family = "quasipoisson")
model.os3 <- glm(num_fake_shares_os ~ party + age + female + black + educ + faminc + num_posts, 
                 weights = weight_svypop_w3, data = dat, family = "quasipoisson")
model.os4 <- glm(num_fake_shares_os ~ party + age + female + black + educ + faminc, 
                 weights = weight_svypop_w3, data = dat, family = "quasipoisson")

varnames <- c("(Intercept)", "Very Liberal", "Liberal",
              "Moderate", "Conservative", "Very Conservative", "Democrat",  "Republican", "Independent", "Other", 
              "Age: 30-44", "Age: 45-65", "Age: Over 65", "Female", "Black", 
              "Education", "Income", "Number of FB posts") 

stargazer(model.os1, model.os2, model.os3, model.os4,
          style = "apsr", align = TRUE, 
          title = "Determinants of Fake News Sharing on Facebook",
          covariate.labels = c(varnames[-1]),
          dep.var.labels = c("Number of Stories Shared"), 
          omit.stat = c("f", "ser", "rsq"), notes = "Quasi-Poisson models.")


## TABLE S13

model.cp1 <- glm(num_fake_shares ~ ideology + age + female + black + educ + faminc + knowl + num_posts, 
                 weights = weight_svypop_w3, data = dat, family = "quasipoisson")
model.cp2 <- glm(num_fake_shares_ag ~ ideology + age + female + black + educ + faminc + knowl + num_posts, 
                 weights = weight_svypop_w3, data = dat, family = "quasipoisson")
model.cp1a <- glm(num_fake_shares ~ ideology + age + female + black + educ + faminc + knowl, 
                  weights = weight_svypop_w3, data = dat, family = "quasipoisson")
model.cp2a <- glm(num_fake_shares_ag ~ ideology + age + female + black + educ + faminc + knowl, 
                  weights = weight_svypop_w3, data = dat, family = "quasipoisson")

varnames <- c("(Intercept)", "Very Liberal", "Liberal",
              "Moderate", "Conservative", "Very Conservative", 
              "Age: 30-44", "Age: 45-65", "Age: Over 65", "Female", "Black", 
              "Education", "Income", "Knowledge", "Number of FB posts") 

stargazer(model.cp1a, model.cp1, model.cp2a, model.cp2,
          style = "apsr", align = TRUE, 
          title = "Determinants of Fake News Sharing on Facebook",
          covariate.labels = c(varnames[-1]),
          dep.var.labels = c("Number of Stories Shared", "Number of Stories Shared (A\\&G)"), 
          keep.stat = c("n", "null.dev", "res.dev"), 
          notes = "Quasi-Poisson models.")

## TABLE S14

model.cp1 <- glm(num_hn_shares ~ ideology + age + female + black + educ + faminc + num_posts, 
                 weights = weight_svypop_w3, data = dat, family = "quasipoisson")
model.cp2 <- glm(num_hn_shares ~ party + age + female + black + educ + faminc + num_posts, 
                 weights = weight_svypop_w3, data = dat, family = "quasipoisson")

varnames <- c("Very Liberal", "Liberal",
              "Moderate", "Conservative", "Very Conservative", "Democrat", "Republican", "Independent", "Other",
              "Age: 30-44", "Age: 45-65", "Age: Over 65", "Female", "Black", 
              "Education", "Income", "Number of links posted", "Constant")

stargazer(model.cp1, model.cp2, 
          style = "apsr", align = TRUE, 
          title = "Determinants of Hard News Sharing on Facebook",
          covariate.labels = c(varnames),
          dep.var.labels = "Number of Stories Shared", 
          keep.stat = c("n", "null.dev", "res.dev"), 
          notes = "Quasi-Poisson models.")

## FIGURE S1

ggdat <- dat %>% filter(!is.na(age)) %>% group_by(age) %>% summarize(mean_fakes = mean(num_fake_shares_ag_url, na.rm = TRUE),
                                                                     mean_fakes_w = weighted.mean(num_fake_shares_ag_url, w = weight_svypop_w3, na.rm = TRUE),
                                                                     sd_fakes = sd(num_fake_shares_ag_url, na.rm = TRUE)/sqrt(n()),
                                                                     sd_fakes_w = sqrt(weighted.var(num_fake_shares_ag_url, w = weight_svypop_w3, na.rm = TRUE)/n()))

p <- ggplot(filter(ggdat, !is.na(age)), aes(x = age, y = mean_fakes_w)) + geom_bar(stat = "identity") +
  geom_errorbar(aes(ymin = mean_fakes_w-1.96*sd_fakes_w, ymax = mean_fakes_w+1.96*sd_fakes_w), width = .1, col = "darkgray")
pq <- p + ggtitle("") + ylab("Mean number of fake news shares") + 
  xlab("Age") + theme(legend.position = "bottom", legend.title = element_blank()) #+

pq + theme_ipsum_rc()
